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Abstract 

An extension of a recently-developed linear thermoelastic theory for multiphase periodic materials 
is presented which admits inelastic behavior of the constituent phases. The extended theory is capable 
of accurately estimating both the effective inelastic response of a periodic multiphase composite and 
the local stress and strain fields in the individual phases. The model is presently limited to materials 
characterized by constituent phases that are continuous in one direction, but arbitrarily distributed within 
the repeating unit cell which characterizes the material’s periodic microstmcture. The model’s analytical 
framework is based on the homogenization technique for periodic media, but the method of solution for 
the local displacement and stress fields borrows concepts previously employed by the authors in con- 
structing the higher-order theory for functionally graded materials, in contrast with the standard finite- 
element solution method typically used in conjunction with the homogenization technique. The present 
approach produces a closed-form macroscopic constitutive equation for a periodic multiphase material 
valid for both uniaxial and multiaxial loading. The model’s predictive accuracy in generating both the 
effective inelastic stress-strain response and the local stress and inelastic strain fields is demonstrated 
by comparison with the results of an analytical inelastic solution for the axisymmetric and axial shear 
response of a unidirectional composite based on the concentric cylinder model, and with finite-element 
results for transverse loading. 


1 Introduction 

Micromechanical modeling of multiphase materials with inelastic phases continues to be a challenging 
problem due to the path-dependence of the local field variables which govern the overall or macroscopic 
behavior of the material. The relationship between the local and global variables through the use of Hill’s 
strain or stress concentration tensors, Hill (1963), which define the instantaneous macroscopic response, 
must be established at each instant along the loading path. This is in contrast with the elastic problem which 
requires determination of Hill’s concentration tensors just once. 
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Many different approaches with varying degrees of complexity and limitations have been proposed for 
modeling the inelastic response of multiphase materials, cf. Dvorak (2000) for a recent review article. Mod- 
els for which exact analytical inelastic solutions of the governing field equations are available are limited to 
the concentric cylinder geometry under axisymmetric or axial shear (but not transverse) loading, Pindera et 
al. (1993) and Williams and Pindera (1997). Finite-element analyses of repeating unit cells representative 
of materials with periodic microstructures produce very accurate estimates of the local field variables, albeit 
at a substantial computational cost. The construction of inelastic macroscopic constitutive equations for 
arbitrary loading based on such analyses is not straightforward, making it difficult to embed finite-element 
based models into more general structural analysis procedures. This has given rise to the development of 
approximate models which employ simplifying assumptions on the form of displacement or stress fields 
within the individual phases. The effect of the inelastic response of the matrix phase is often accounted for 
through the average values of the inelastic strain field either in the entire matrix phase or in the subvolumes 
into which the matrix phase is partitioned. 

A number of approximate analyses of a repeating unit cell, which employ subvolume discretization of 
the unit cell to mimic the material’s microstructure, has been developed to deal with the inelastic response 
of periodic multiphase materials. The Generalized Method of Cells, Paley and Aboudi (1992), employs 
a first-order representation of the displacement field in each subcell of the repeating unit cell, producing 
piece-wise uniform strain and stress fields throughout the cell. The method is a generalization of the 
original Method of Cells developed by Aboudi (1982) which is based on a limited cell discretization. In the 
context of unidirectional composites, the generalized version allows accurate and efficient analysis of the 
impact of fiber shape and arrangement on the composite’s inelastic macroscopic response, as demonstrated 
by Arnold et al. (1996). The predictive capability of the method in various applications has been recently 
summarized by Aboudi (1996). However, despite the method’s demonstrated accuracy in modeling the 
inelastic macroscopic response of periodic composites, the accuracy with which local stress and strain fields 
are captured (although acceptable in many applications) is not as good. This, for instance, requires the 
incorporation of additional assumptions and modifications into the model’s framework in order to trace 
the evolution of damage at the local level due, for instance, to fiber breakage and fiber/matrix debonding, 
Bednarcyk and Arnold (2001a, b) and Pahr and Arnold (2001). 

The assumption of piece-wise uniform strain fields has also been employed by Dvorak (1992) in the 
context of a procedure called the Transformation Field Analysis. This approach has recently been gener- 
alized by Chaboche et al. (2001) and demonstrated to capture the local stress and (inelastic) strain fields 
with good accuracy in comparison to finite-element simulations. Methods which employ more accurate field 
representations within the repeating unit cell of a periodic composite, albeit at a significantly enhanced com- 
putational cost, have been developed by Walker et al. (1994) and Fotiu and Nemat-Nasser (1996). These 
are based on Fourier series approximations of the stress and strain fields within the repeating unit cell. In 
practice, the solution of the local field quantities is obtained by discretizing the cell into square or triangular 
subvolumes in which the eigenstrains are assumed to be uniform. 

Most recently, a new method for the analysis of the thermoelastic response of multiphase periodic mate- 
rials, characterized by repeating unit cells with arbitrary microstructures in the plane normal to the direction 
of the continuous reinforcement (admitting fully anisotropic behavior in this plane), was developed in which 
the displacement field within each subvolume of the repeating unit cell was approximated by quadratic func- 
tions expressed in local coordinates, Aboudi et al. (2001). This, in turn, produces linear strain and stress 
fields at the local subvolume level, in contrast with the piece-wise uniform fields employed in the Gener- 
alized Method of Cells and Transformation Field Analysis. The method’s analytical framework is based on 
the homogenization technique (cf. Sanchez-Palencia, 1980; Suquet, 1987; Parton and Kudryavtsev, 1993; 
Kalamkarov and Kolpakov, 1997), but the solution procedure for the local displacement, strain and stress 
fields within each subvolume of the repeating unit cell borrows concepts previously employed by the au- 
thors in constructing the higher-order theory for functionally graded materials, Aboudi et al. (1999). The 
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higher-order subvolume displacement representation relative to that employed in the Generalized Method 
of Cells provides the necessary shear-coupling between the local normal and shear deformation fields and 
the macroscopically applied average strains. This feature, in turn, produces excellent estimates of both the 
effective thermoelastic moduli of periodic composites and the local stress and strain fields in the individual 
constituents, as demonstrated in Aboudi et al. (2001) through comparison with the finite-element results of 
Sun and Vaidya (1996) for the effective moduli and an analytical solution for the local stress fields. There- 
fore, in view of the above features and the same volume discretization methodology as that employed in the 
Generalized Method of Cells, we refer to this new method as High-Fidelity Generalized Method of Cells or 
HFGMC. 

Herein, an extension of this method is presented which accounts for the inelastic behavior of the in- 
dividual phases. The method’s fully analytical nature results in a closed-form expression for the effective 
response of a periodic multiphase material under multiaxial loading given in the form of a macroscopic 
thermoinelastic stress-strain relationship, which can easily be incorporated into a structural analysis code as 
a subroutine. The ease of the repeating unit cell’s construction and the demonstrated accuracy of the theory 
in predicting the average and the local stress and plastic strain fields, facilitates quick, efficient and reliable 
analysis of the impact of a multiphase material’s microstructure on the average and local response. The the- 
ory’s accuracy is validated based on comparison with the analytical concentric cylinder model predictions 
for axisymmetric and axial shear loading, and finite-element predictions for transverse loading. Further, 
comparison with the predictions of the original Generalized Method of Cells is presented in order to illus- 
trate the advantages derived from the new approach. This comparison also provides an explanation for the 
good accuracy of the Generalized Method of Cells at the macroscopic level despite poor estimates of certain 
stress components at the microscopic level. 

2 Theoretical Framework 

Consider a multiphase composite wherein the microstructure is periodically distributed in the plane x 2 — 
x% defined by the global coordinates (^2,^3), see Fig. 1 where the repeating unit cell used to construct 
the periodic array is highlighted. In the framework of the homogenization method, the displacements are 
asymptotically expanded as follows 


Ui (x, y ) = u Qi (x, y ) + 8 u u (x, y ) + 6 2 u 2 i (x, y ) + . . . 


( 1 ) 


where x = (x\, x<i, #3) are the macroscopic (global) coordinates, and y = (2/1, 2/2? 2/3) are the microscopic 
(local) coordinates that are defined with respect to the repeating unit cell. The material’s periodicity imposes 
the constraint u a ifa y) = ^(x, y+n p d p ) on the different-order terms u a i (a = 0, 1, 2, ...) in Eq. (1), 
where n p are arbitrary integer numbers and the constant vectors d p characterize the material’s periodicity. 
In addition, the size of the unit cell is further assumed to be much smaller than the size of the body so 
that the relation between the global and local systems is y\ = Xi/6, where 6 is a small scaling parameter 
characterizing the size of the unit cell. This implies that a movement of order unity on the local scale 
corresponds to a very small movement on the global scale. 

Employing the following relation in evaluating the derivative of a field quantity: 


d did 

dxi dxi ^ 8 dyi 


( 2 ) 
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Figure 1. A multiphase composite with a periodic microstructure in the x 2 — £3 plane characterized by a 
repeating unit cell (highlighted). 


the strain components are determined from the displacement expansion ( 1 ) in the following form 

eij = €ij(x) + €ij(x, y) + 0(6) (3) 

where the average and fluctuating strains, i,j (x) and (x, y ), are given by global and local strain-displacement 

relations as follows 


13 2 dxj dxi * 


~ , , 1 f dui duj 


( 4 ) 


In the above equation, U{ are displacement components in the homogenized region (at the continuum scale) 
and hence are not functions of the local coordinates y, and uu = U{ are the fluctuating displacements. It 
can be easily shown that 


vJ'M-vJc 


^ij + ^ij ) dVy — 6 


V 


where V y is the volume of the repeating unit cell. This follows directly from the periodicity of the fluctuating 
strain, implying that the average of the fluctuating strain taken over the unit repeating cell vanishes. 

Using Eq. (3), one can readily represent the displacements in the form 


Ui(x, y) = eijXj +Ui + 0(6 2 ) (5) 

where e^xj represents the contribution of the average (homogenized) strain to the total displacement field 
and Ui represents the fluctuating displacement field. This representation will be employed in constructing 
an approximate displacement field for the solution of the cell problem discussed in Section 2.2. 

For given values of the average strains e^, the unknown fluctuating displacements are governed by the 
equilibrium equations subject to periodic boundary conditions imposed on the displacement and traction 
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components that are prescribed at the boundaries of the repeating unit cell. In addition to these boundary 
conditions one needs to impose the continuity of displacements and tractions at the internal interfaces be- 
tween the phases that fill the repeating unit cell. The manner of solving the governing field equations for the 
fluctuating displacements in the repeating unit cell based on the representation given by Eq. (5) is described 
next. 


2.1 Definition of the cell problem 


We perfom the analysis on the repeating unit cell which occupies the region 0 < 2/2 < 77, 0 < 2/3 < £ 
specified in terms of the local coordinates (2/2 , 2/3) • The microstructure of the repeating unit cell in the 2/2 — 2/3 
plane in the given region is discretized into N q and N r internal or generic cells. Figure 2(a) illustrates how 
the repeating unit cell highlighted in Fig. 1 could be discretized. In addition, every generic cell consists 
of four subcells designated by the pair (^ 7 ) where each index takes the value 1 or 2 which indicates the 
relative position of the given subcell along the 2/2 and 2/3 axis, respectively, see Fig. 2(b). The indices 
q and r, whose ranges are q = 1, 2, ..., N q and r = 1,2, ..., N r , identify the generic cell in the 2/2 — 2/3 
plane. The dimensions of the generic cell along the 2/2 and 2/3 axes are h^\h^ and l^\ \ such that 

H = YpliPP + ^ = i(P + 4 0 ). This manner of discretizing a periodic material’s 

microstructure has also been employed in constructing the original Generalized Method of Cells (Paley and 
Aboudi, 1992). The construction of the higher-order theory for functionally graded materials, characterized 
by spatially variable microstructures without a definable repeating unit cell, is also based on such volume 
discretization (Aboudi et al. 1999). 

Given an applied macroscopic loading specified by the average strains e^-, an approximate solution for 
the displacement field within each (^ 7 ) subcell of the N q N r generic cells is constructed based on volumetric 
averaging of the equilibrium equations together with the imposition of periodic boundary conditions, and 
both the displacement and traction continuity conditions, in an average sense between the cells and subcells 
used to characterize the material’s microstructure. The equilibrium equations for the (/Ty) subcell occupying 
the region | yp \< / 2 , | yP |< p / 2 that the displacement field must satisfy are given by 


+ d 3 ag 7) = 0 je 1,2,3 ( 6 ) 

where O 2 = O/dy^' and 83 = d/d y^\ The subcell stress components are related to the subcell strains 
through generalized Hooke’s Faw which includes spatially uniform thermal loading characterized by the 
temperature deviation AT from a reference temperature and inelastic effects, 



a, 


(07) 




(07) 
kl 


rnP7; (API) . 

0 ijkl V e kl e 


/(07) 

L kl 


— 6 


T(0 7 )> 
kl > 


(7) 


where Cpp are the elements of the stiffness tensor of the (/Ty) subcell , e^p are the total strains, and , 

are the inelastic and thermal strains in the subcell, with no summation implied by repeated Greek 
letters in the above and henceforth. The inelastic strains are derived from the chosen inelastic constitutive 
model. In this paper, we consider either elastic orthotropic materials (characterized by nine independent 
Cijki elements) or inelastic materials which are isotropic in both elastic and inelastic domains. Hence, Eq. 
(7) reduces to (assuming incompressibility of the inelastic strains) 


A07) _ ^(07) .(07) _ 9. (07) 7(07) _ T(J3j) 

a ij ° ijkl e kl Hi a ii 






( 8 ) 


where is the elastic shear modulus of the material filling the given subcell (^ 7 ), and the term 
henceforth referred to as thermal stress, stands for the thermal contribution 


a, 


7X07) _ r (07) 




rlP7U^ 

ij 


(9) 
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where T^ 7 ^ are the thermal stress coefficients. We relate the subcell displacement field, which forms the 
basis for the solution of the equlibrium equations, to the subcell stresses through the Hooke’s Law (7) in 
conjunction with the subcell strain-displacement equations. Taking into account Eqs. (3) and (4), the total 
strains in the subcell (/Jy) are given by 


(£7) 

ij 


= cy + 2 ( diU : 


m + djuf 1 ^ 


( 10 ) 


where 8 \ = 0 and 82 , 83 have been defined previously. 


A 


y 2 


L 



4 (r) e 2 (r) 


generic cell (q,r) 


(b) 


Figure 2. (a) Volume discretization of the repeating unit cell employed in the present model, (b) generic cell 
within the repeating unit cell. 

2.2 Method of solution for the cell problem 

We begin the solution of the equilibrium equations by approximating the fluctuating displacements in each 
subcell using a quadratic expansion in terms of local coordinates {y^\ centered at the subcell’s mid- 
point. This is in sharp contrast with the original Generalized Method of Cells where the employed displace- 
ment expansion was linear as a result of which the coupling between the local normal and shear effects was 
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lost. A higher-order representation of the fluctuating displacement field is necessary in order to capture the 
local effects created by the field gradients and the microstructure of the composite. This was demonstrated 
within the framework of the elastic version of the High-Fidelity Generalized Method of Cells in the context 
of the problem of an infinite plate with a circular inclusion subjected to a far-field uniform stress (Aboudi et 
al. 2001). Both qualitative and quantitative aspects of the local stress field within the inclusion and its vicin- 
ity, as well as the asymptotic character of the stress field away from the inclusion, were captured with very 
good accuracy upon comparison with the exact analytical solution. Such local effects cannot be captured by 
the original Generalized Method of Cells due the first-order (linear) displacement representation within each 
subcell (which produces piece-wise uniform strain and stress fields within the repeating unit cell), and the 
imposition of traction continuity conditions at the subcell interfaces in an average sense. As a consequence, 
tractions along each column of subcells spanning the representative volume element in the directions as- 
sociated with the particular traction components are uniform, precluding the posssibility of modeling local 
effects due to matrix-inclusion as well as adjacent inclusion interactions (Pindera and Bednarcyk, 1999). 
This shortcoming has been eliminated by the present theory as previously demonstrated in the elastic case 
and as will be demonstrated in the inelastic case considered in the current contribution. 

Following the general displacement representation for periodic media given by Eq. (5), the subcell 
displacement approximation in the High-Fidelity Generalized Method of Cells has the form (omitting the 
cell label (q, r), see Fig. 2b), 


7,(^) — . X . _L W^) _L yW W^) 

U l e lj X J ' ^(OO) + y 2 VV 1(10) 


1(10) +2/3 *^1(01) 


+ 




fc(«) 2 Ar) 2 

;<3 Sf 2 -V)'«7( 3s 3 7l2 -V 


4 ' ”1(02) 


+ 


t2] x ] T ''2(00) + y 2 ''2(10) + W3 ''2(01) 

h (^) 2 1 /(r)2 


- 




> '' 2 ( 02 ) 


( 11 ) 


( 12 ) 


u. 


(37) 


where W^oo J are the fluctuating volume-averaged displacements, and ( 7 = 1, 2, 3) are the higher- 

order terms. The number of these unknown microvariables that describe the fluctuating displacements in 
the cell (#, r) is 60. These microvariables are determined by satisfying the equilibrium equations (6) in 
a volumetric sense, the interfacial continuity conditions (both displacements and tractions) on the faces 
separating adjacent subcells and cells in an integral sense, 


6 3 J X J I ' ' 3(00) + y 2 ^(lCT) ' y 3 VVr - 


+(10) 


.-7(7)1 


(^7) 

3(01) 


+ 


^(3y? )2 - 


_3 

A ' ' ' 3(20) 


+ i(3pf 2 - 


(r)2 


r(3 7) 


' ''3(02) 


(13) 


h/ 2 


h/ 2 


Z 7 /2 


JM) |(9- r ) AJn) 
^ VV= hl /2 dy z 


1 7/2 


(2 7 ) ,(g,r) , h) 

i V?=-h 2 /2 dy z 


hp/2 

J 

-h, 3/2 


r (/31) I (q,r) 


/2^ 


( 3 ) 


hp/2 


~h, 3/2 


,(32) i(9,r) 


'3 j 


4 2) =-V 2 


dy f 
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Z 7 /2 Z 7 / 2 

f (17) I (9+1.0 ^-( 7 ) _ f (27) |(9.0 jij(7) 

J Lw=- fcl /2 dy 3 - J °2j L( 2 ) =fe2/2 ^3 


J ^=- hl/ 2 Uy * J 'y?=h 2 / 2 Uy * 

-l-y/2 -ly/ 2 

hp/ 2 ft/g/2 

/ 4fi|0> i/2 #f = / 

ft /3 /2 ft /3 j2 

ly/2 ly/2 

f rr(l7) |( g - r ) WyM - / 

./ * W^i /2 ^ - J 


a m \ {q ’ r) dv W) 
^ l y^=l 2 /2 aV2 


0/2 

ft/ 3/2 

I 

~hp / 2 
0/2 


/ rr(27)|(9,0 ^(7) 

J U i \y?=- h 2 / 2 dy * 

' 0/2 


^/2 

~(/31) i(9,0 j^(0) _ f 57 (^ 2 ) i(9,0 ^(/3) 

W^i/2^ 2 ~ J Ui V 3 2 ) =-h / 2 ^ 

-ft/ 3/2 


0/2 

[ 57(17) 1(9+00 ^(7) 


Z 7 /2 




fytf / 2 


I 

-ft/ 3/2 


u. 


(P !) |(9,r+l) 
0 ( 1 ) 


2/3 —— 0/2 




ft/3/2 

■s 

-ft/ 3/2 


</32) 1 ( 9,0 


^ 2 ) = 0/2 




(14) 


and the periodic boundary conditions at 1/2 = 0,H and 2/3 = 0 , L expressed in the local subcell coordinates 
for the appropriate subcells within the boundary generic cells 


fc 7/ 

/ 


- 0/2 


.(I 7 ) |( 1,0 


2 / 


4 1} =' 


fti/2 


#3 7) 



(7 


(27) |OVO 

2 / l ^ 2 ) =ft 2 /2 



ft/3/2 

/ 

-ft/3/2 

Z 7 /2 


ft/3/2 


^ l y « = _ h/2 ^2 


J 31 »?-V2 ” 2 


u. 


(17) 1(1,0 

i L(l) _ 


- 0/2 

ft/3/2 

/ 

-ft/3/2 


1/2 =— fti/2 




( 7 ) 


« (/31) l (9,1) du (/3) 

l yW = _ Zl/2 ^2 


-ft/ 3/2 

ly/2 

-I'' 

~ 0/2 

ft/ 3/2 

-S 

-ft/ 3/2 


(27) I (Nq,r) , h) 


i U2)_ 


^3 


vr=h 2 / 2 


~(/32) I (q,N r ) Uf3) 
1 '^ 2) = 0/2 2/2 


(15) 


where j = 1, 2, 3. These periodicity conditions ensure that the repeating unit cell is an intergral part of 
the periodic array representing the entire multiphase material, thereby explicity accounting for the interac- 
tion with the adjacent unit cells. As shown in the following development, application of these conditions 
produces the correct number of equations for the determination of the unknown microvariables. 
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In the perfectly elastic case, the above displacement expansions produce linear variations in strains and 
stresses at each point within the subcell. In the presence of inelastic effects, however, a linear strain gener- 
ated by these equations does not imply the linearity of the stress field due to the path-dependent deformation. 
Thus the displacement field microvariables must depend implicitly on the inelastic strain distributions, giv- 
ing rise to a higher-order stress field than the linear strain field generated from the assumed displacement 
field representation. In the presence of inelastic effects, this higher-order stress field is represented by a 
higher-order Legendre polynomial expansion in the local coordinates. Therefore, the strain field gener- 
ated from the assumed displacement field, and the resulting stress field, must also be expressed in terms of 
Legendre polynomials: 


f (Pi) 

ij 


^1 + 2 m) (1 + 2 n)4j(l, n ) J 


( 16 ) 


“I - 2?7l)(l 27?/) 


P m (Cf)Pn(d 7) ) 


where the non-dimensional variables defined in the interval — 1 < < 1, are given in terms of the 

local subcell coordinates as (P = yP /(h^ /2), and (P = yp /{lP /2). For the given displacement 
field representation, the upper limits on the summations in Eq. (16) become 1. The upper limits on the 
summations in Eq. (17) are chosen so that an accurate representation of the stress field (which depends on 
the amount of the inelastic flow) is obtained within each subcell, Aboudi et al. (1999). The coefficients 
e>! d p x, t^P x in the above expansions are determined as described below. 

The strain coefficients pjp n ^ in subcell (/3j) of cell (#, r) are explicitly determined in terms of the 
displacement field using the orthogonal properties of Legendre polynomials. The non-zero components are 
given as follows (omitting (q,r)) 

e ( ^ } -en 

e ii(o,o) — en 

p (P 7) p 1 TX/iPl) 

e 22(0,0) 622 * / ^2(10) 

p (^7 ) _ u(<l) rxr{Pl) 

e 22(l,0) 2 l P ^2(20) 

g 1 ^7) _ g 1 w(07) 

e 33(0,0) _ e33 ^ ^3(01) 

{Pi) — /( r )w(/37) 

e 33(0,l) 2 7 3(02) 

£ {Pl) _ g I L W Vh) I w(^7) \ 

® 23(0,0) — 623 + 2 2 ( 01 ) + ^3(10)' 


p {Pi) _ y^Mu/iPi) 

e 23(l,0) “ 4 ,1 P 1 1 3(20) 


-23(0,1) 


= —Wwtt 


piPl) _ p.„ 1 

e i3(0,0) 613 + 2 ^1(01) 


"13(0,1) 


= 
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p (Pi) ~p _i_ 

e i2(0,0) — e 12 + 2 ^1(10) 

p (Pi) — 

e i2(l,0) 4 hp ^1(20) (1°) 

It should be noted that 4f(c? 0 ) P rov ide the average strains in subcell (f3j) of cell (q, r). 

The stress coefficients in subcell (f3j) of cell (q, r ) are expressed in terms of the strain co- 

efficients, thermal stress and tne unknown inelastic strain distributions, by first substituting the Legendre 
polynomial representations for the total strain and stress into the constitutive equations, Eq. (8), and then 
utilizing the orthogonality of Legendre polynomials: 


r {Pi) _ -n[Pi) e ^) _ 

' ij(m,n) ^ ijkl lh ij(m,n) 


CrJj^SmoSon 


The R 


!j(m n) terms represent the inelastic stress distributions calculated in the following manner 


lj(m,n) 


= V^V(2m + l)(2n + l) JJ^PM^Pn^)^^, 




Note that in both Eqs. (19) and (20) the cell labeling (#, r) has been omitted. 

In the course of satisfying the equilibrium equations in a volumetric sense, it is convenient to define the 
following stress quantities: 


r qiPl) 1 (q,r) 1 

L ij(m,n)i 7 (<?) 7 (r) 

ri ^ Lsy 


r hf/2 r iP/2 


I -h™/ 2 J-l^/2 ° r * 3 




(.y^Tdyfdy? 


For m = n = 0, Eq. (21) provides the average stresses in the subcell, whereas for other values of (m, n) 
higher-order stresses are obtained that are needed to describe the governing field equations of the contin- 
uum. These stress quantities can be evaluated explicitly in terms of the unknown coefficients by 

performing the required volume integration upon substituting Eqs. (8), (10) and (1 1)-(13) in Eq. (21). This 
yields the following non- vanishing zeroth-order and first-order stress components in terms of the unknown 
coefficients in the displacement field expansion (omitting (<7, r)): 


APi) _ 


"11(0,0) 


= C^en+C- 


(&y)(uAP'y) 


+ C22) + C'i3 7) (W / 3foi') 


e 33 )-rf 7) 


AT — 

L ^" L -”-11(0,0) 


c(P'y) \ ,(<?) P037) 

°11(1,0) 4 n (3 °12 VV 2(2Q) 2v / 3 ,t/? rtll ( 1 ’°) 

aWl) _ I/(r)2^(/37)ii/(^7) _ _J_/(r) r>(J3y) 
O ll(0,l) 4*7 °13 ^3(02) 2^/3 7 rlll ( 0 ’ 1 ) 

with similar expressions for the other normal stress components, and 


7 ) (Op 1 wW 1 txt (^7) ■ 

*23(0,0) ^44 V^ e 23 + rr 2 ( 01 ) + ^(lop 


*-23(0,0) 


= -h 


(q) 2 r<(Py)uAPi) 


'23(1,0) 4 73 


1 u(q) r>(Jh) 

,— ILn -TLoo/t , 


3(20) 2 y/3 P “23d, 0) 


c(/?7) _ ^ / (r)2/~<{Pl) ’uriPl) _ \ Ar) r>(Pl) 

°23(0,1) 4 ( '7 °44 ^2(02) “"l 0 ’ 1 ) 


? 13(0,0) 


= Qs (2e 13 + W"uofO - R 


13(0,0) 


NASA/TM— 2002-21 1469 


10 



q(97) _ 3 7 (r)2 ^(97)T;r/(97) 3 ;(r)p(97) 

°13(0,1) _ 4 7 U 55 ”1(02) 2^3 7 ^(O’ 1 ) 

(29) 

c(97) _ p>(97)/ 9? _i_w(97)\ d(97) 

°12(0,0) — u 66 ( Z£ 12 “I" 11 1(10)) 2l 12(0,0) 

(30) 

c( 97) _ ^ , (g)2^,(^7)Tr r (/37) 1 .(«) d(97) 

°12(1,0) 4'^ °66 ”1(20) 2 Vo ^ rll2 ( 1 ’ 0 ) 

(31) 


where contracted notation has been employed for the stiffness elements . 

Subsequently, satisfaction of the zeroth, first, and second moments of the equilibrium equations (6) re- 
sults in the following 12 relations among the the volume-averaged first-order stresses n ^ in the different 
subcells (/7y) of the (q, r ) cell, after lengthy algebraic manipulations 

[Sg£o)/*S + S mli)/W q ’ r) = 0 j = !. 2. 3 (32) 

The continuity of tractions at the subcell interfaces within each interior cell, and between adjacent cells, 
imposed in an average sense, can be shown to reduce to the following relations which are obtained from the 
first four equations of (14) after some algebraic manipulations, 


[ 12 ^2i( 1,0) 


+ 5, 


(27) 

2i(0,0) 


- <5Sg(i ) ,0)/fe] tar| - [sg(io) + <3 S w)/'* 2] <,_1 -’' ) = ° 


(33) 


f_c(h) I 1 c( 2 7 ) _ 00 ( 27 ) /L ](g,r) , lrq( 2 7) 

L O 2j(0,0) 2 2j(0,0) OO 2j(l,0)/ A ^2J ■+■ 2 L°2 j(0,0) 


+ 65, 


+ s vm - <(0,i)/fa] <4r| - I s S,0) + 6 -sS, i)/y = 0 


492) 


?(/») 


092) 


(27) 

2j(1.0) 

?092) 


/fc 2 ] (9 “ 1,r) = 0 




1 ^ o(9 2) 

3f (0,0) ^ 2 3j(0,0) 


(91) 


■35. 


(92) n i(9,r) , I r o(92) 
3i(0,l)/ t2 J ^ 2 ^ 


1 « c(9 2 ) 

3f (0,0) + D,5 3i(0,l) 


/Z 2 ] (<? ’ r - 1) = 0 


(34) 

(35) 

(36) 


where j = 1, 2, and 3. The details of derivation of equations similar to Eqs. (32)-(36) have been provided 
by Aboudi et al. (1996) for a more general case of a microstructure containing periodic inclusions in the 
out-of-plane direction. 

Equations (33)-(36) provide 24 additional relations among the zeroth-order and first-order stresses. 
These relations together with Eq. (32), can be expressed in terms of the unknown coefficients W^\ 
by making use of Eqs. (22)-(31), providing a total of 36 of the required 60 equations necessary for the 
determination of these coefficients in the cell (g, r). 

The additional 24 relations necessary to determine the unknown coefficients in the displacement field 
expansion are subsequently obtained by imposing displacement continuity conditions on an average basis at 
each subcell and cell interface. This produces, upon use of the last four equations of (14), 


K‘4 + 


*(io) 


+ -J 


= [wf4 - 7*44 + 7<4i <w ' ) 


i(00) 


*( 10 ) 


i(20)J 


(37) 


<4 + + \hi <4,] h, ’' ): = I<4 - + ^f<4)] <,+1 -’' ) 
[<4 + = <4 - + 7727’ 

[74 + + \‘i<<Sr y - <4 - + 70 ( * r+1) 


(38) 

(39) 

(40) 


where i = 1, 2, and 3, which comprise the required additional 24 relations. 

The equilibrium relations, Eqs. (32), together with the traction and displacement continuity conditions 
Eqs. (33)-(36) and (37)-(40), respectively, form 60 equations in the 60 unknowns W^n) which govern the 
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equilibrium of a subcell (^ 7 ) within an interior cell (< 7 , r); q = 2, ..., N q — 1, r = 2, 7V r — 1. For the 
boundary cells q = 1 , and r = 1, iV r a different treatment must be applied. 

For the boundary cell (1, r), the above relations are operative, except Eqs. (33) and (34), which follow 
from the continuity of tractions between a given cell and the preceding one (since cell ( 0 , r) does not exist). 
These 12 equations must be replaced by the conditions of continuity of tractions at the interior interfaces 
within cell ( 1 , r) (imposed in the average sense), and by the conditions that the fluctuating displacements 
are periodic. It follows, by employing the third set of relations in (15), that 


KS, - |ai<4 + = [<?> + wtf> + 


(41) 


where i — 1,2, and 3. Both conditions provide the required 12 relations to be used for cell (1, r). 

For the boundary cell (N q , r), the previously derived governing equations are operative except for the 
6 relations given by Eqs. (38), which are obviously not applicable (since cell ( N q + 1, r) does not exist). 
These are replaced by the conditions that the tractions are periodic. Thus the first set of relations in (15), 
reproduced below for convenience, provides the 6 equations to be used for cell (N q , r ) 


h! 2 


ii! 2 


cr. 


(dv) 1 ( 1,0 


2 j 


V2 )= - 


M/2 


#3 7) 


111 2 


111 2 


a 


(27) I (N q ,r) 

2j v { ?=h 2 l 2 



(42) 


where the stresses are given by Eq. ( 8 ). 

Similar treatments hold for boundary cells (#, 1) and (q,N r ). Thus the 12 equations (35)-(36) are 
obviously not applicable in cell (< 7 , 1 ) and should be replaced by the conditions of continuity of tractions 
at the interior interfaces within this cell (imposed in the average sense), and by the conditions that the 
fluctuating displacements are periodic. The latter yield according to the fourth set of relations in (15) 


\w m 

r v i( 00 ) 


\hwr 


<( 01 ) 


+ 



(43) 


For the boundary cell (q, Ay), Eqs. (40) are not operative and they should be replaced by the periodicity of 
tractions which is given by the second set of relations in (15), reproduced below for convenience, 


h, 3/2 

f _031) i(9.1) 

J l i4 1) =- 

-h, 3/2 


hp/ 2 


h/2 


dy { 2 0) 


-1 

-h, 3/2 


.032) Aq,N r ) 
'yW=l 2/ 2 


dy { 2 0) 


(44) 


Consequently, the governing equations for the interior and boundary cells form a system of 60A^A^ r 
algebraic equations in the unknown coefficients . The final form of this system of equations can be 

symbolically represented by 

KU = f + g (45) 

where the structural stiffness matrix K contains information on the geometry and thermomechanical prop- 
erties of the materials within the individual subcells (^ 7 ) of the cells comprising the multiphase periodic 
composite. The displacement vector U contains the unknown displacement coefficients in each subcell, i.e., 

U =[<>,.. (46) 

where in subcell (^ 7 ) of cell (< 7 , r) these coefficients are 

ug 7) = [w i{ 00 ), Wi(io), W m) , w im , w im }^ * = 1, 2, 3 
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The mechanical force vector f contains information on the applied average strains and the imposed tem- 
perature deviation AT, and the inelastic force vector g appearing on the right-hand side of Eq. (45) contains 
the inelastic effects given in terms of the integrals of the inelastic stress distributions that are represented by 
the coefficients ^ . These integrals depend implicitly on the elements of the displacement coefficient 
vector U, requiring an incremental procedure for the solution of Eq. (45) at each point along the loading 
path. 

A careful check of the preceding equations reveals that the equations that govern the local normal and 
in-plane (2-3) shear deformations are coupled, thus providing the necessary shear coupling effects. On 
the other hand, these equations are not coupled to the axial shear deformations (1-2 and 1-3). Thus the 
above system, Eq. (45), can be decoupled in practical applications and solved for the normal and transverse 
shear deformations (with AQN q N r algebraic equations) separately from the axial shear deformations (with 
2QN q N r algebraic equations). 


2.3 Global constitutive relations 


Once the solution U for a given set of average strains e has been established, we can determine, in particular, 
the average strains in subcell (/Ty) of the cell (q, r ) given by (18). The average stress components 

su bcell (Pj) of the cell (g,r) are given by Eqs. (22), (25), (28) and (30). They can be 
assembled in a compact form as follows 

[s® ( * r) = [c^egg - R(oo) - T^AT}^ (47) 

The equation relating the average total, plastic and thermal subcell strains, and macroscopically applied 
strains is obtained by generalizing the localization relation given in Aboudi et al . (2001) for elastic phases 
in the following manner 

[ e (0,0)] (9,r) = [A (/37) e + D^ 7 )]^) (48) 

where [A^ 7 )]^) is the mechanical strain concentration matrix of the subcell (/Ty), and [D^T^ r ) is a 
vector that involves current thermo-inelastic effects in the subcell. In the absence of thermal and inelastic 
effects this vector vanishes, and we can readily determine from (48) the mechanical strain concentration 
matrix [A^ 7 )]^ 7 *) by solving the system (45) six consecutive times upon imposing a single non-zero com- 
ponent of e one at a time. 

The thermo-inelastic analysis is performed in conjunction with an incremental procedure by imposing 
a spatially uniform temperature AT and the applied macroscopic strain e in a stepwise manner. Thus for a 
given value of applied thermomechanical loading, the average strains in the subcell are obtained 

from the solution of Eq. (45), and hence the matrix [D^ 7 )]^ ,r ) from (48) at the current loading level. 
Substitution of (48) into (47) yields 

[ S (S (, ’ r) = [C ( ^ 7) (A ( ^ 7) e + D (/?7) ) - rJJJJ - T^AT] {q ’ r) (49) 

The average stress in the multiphase periodic composite is determined from 


(T 


1 

Tl 


Ng Nr 2 

EE E 


q=l r 1 (3, 7=1 


(50) 


Consequently, Eqs. (49)-(50) establish the effective constitutive law of the multiphase thermo-elastic 
composite in the form 

<t = C* 6 — (a 1 + a T ) (51) 
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where C* is the effective elastic sti ffn ess matrix of the composite which is given by 


C* 


TT H hfl^[ 


HL 


q = 1 r = 1 /3,7=1 


( 52 ) 


and a 1 and cr T denote the overall (macroscopic) inelastic and thermal stresses in the composite whose sum 
is given by 


N q N r 2 

»'+* T = ^EE E [C (/?7) D (/?7) - AT] ( «’ r) (53) 

q = 1 r = 1 /3,7=1 

Owing to the implicit dependence of on the inelastic and thermal effects, the right hand side of 
(53) involves the combined inelastic and thermal contributions. It is very convenient to utilize the result 
obtained by Levin (1967) to separate the global inelastic and thermal effects. The global thermal stress in 
the multiphase composite a T = T* AT ( T* is related to the effective coefficients of thermal expansion a* 
of the composite by T* = C*a* ) is given in accordance with the Levin formula in terms of the mechanical 
strain concentration matrices and the thermal stress vector in the individual phases by 


<j = 


AT 

Wl 


N q N r 


<7=1 r = 1 /3, 7=1 


(54) 


where [A tr ^ 7 )]^ ,r ) is the transpose of the mechanical strain concentration matrix [A^ 7 )]^ ,r ) of subcell 
(/? 7 ) within cell (q,r). This provides an additional check on the consistency of the proposed approach. 
Consequently, by utilizing Eq. (54) the overall inelastic stress of the composite can be readily obtained from 
(53), so that it can be represented in a closed-form manner, namely, 


N q N r 2 

a 1 = T T h 0 )l ^ r) - r[(7) + (A tr ^ - I)r^ 7) A T] (9 ’ r ) (55) 

<7=1 r = 1 /3,7=1 

We remark that, in the present framework, we have chosen to use the total formulation approximation 
of the field quantities which, in turn, results in the total form of the effective stress-strain relations for the 
multiphase material given by (51). The actual integration of the inelastic effects represented by the vector 
[D(^ 7 )]( 9,r ) depends on the chosen inelastic constitutive model for the individual phases. In the present 
paper, these effects are integrated using an iterative procedure within an incremental framework described 
in the next section for the incremental plasticity theory representation of the inelastic phase response. In the 
case of rate-dependent phase constitutive models, which the High-Fidelity Generalized Method of Cells’ 
framework can naturally accomodate, either an explicit or an implicit integration scheme may be employed. 
This is in contrast with the approach employed by Paley and Aboudi (1992) where both total and rate 
formulations were employed in approximating the subcell field quantities. Similarly, rate formulations were 
employed by Fotiu and Nemat-Nasser (1996) and Walker et al . (1994). Chaboche et al. (2001) employed 
both rate (or tangential) and secant formulations in their extension of the Transformation Field Analysis 
which also included the presence of evolving damage modeled using the continuum damage mechanics 
approach. The present theoretical framework can also be reformulated in terms of rates, however for stand- 
alone applications there is no clear advantage of the total vs rate approach. 
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3 Numerical Results and Discussion 


The focus of the extended High-Fidelity Generalized Method of Cells’ validation is both the accurate predic- 
tion of the macroscopic inelastic response of a unidirectional metal matrix composite and the evolution of the 
local stress and inelastic strain fields at different macroscopic loading levels. The multiple concentric cylin- 
der model, briefly described in the Appendix for completeness, is employed to validate the High-Fidelity 
Generalized Method of Cells’ predictive capability under axisymmetric loading due to a spatially uniform 
temperature change (Pindera et al ., 1993), as well as under axial shear loading (Williams and Pindera 1997). 
This model is chosen because it provides fully analytical solutions for the local stress and inelastic strain 
fields, which satisfy the external boundary conditions together with the fiber/matrix interfacial displacement 
and traction continuity conditions under such loading. Therefore, these solutions are exact and provide the 
effective elastic and inelastic response of unidirectional composites with arbitrary fiber content. However, 
because the solutions are based on a particular geometry and the concept of a representative volume element 
which employs homogeneous boundary conditions to obtain solutions for the local fields, the functional 
form of the local fields is not affected by the presence of adjacent fibers. At dilute fiber volume fractions, 
the presence of adjacent fibers and the actual fiber array geometry has very little effect on both the local 
field quantities and the macroscopic response, and therefore the High-Fidelity Generalized Method of Cells 
which is based on periodicity (and the use of periodic rather than homogeneous boundary conditions) can 
be compared directly with the multiple concentric cylinder model predictions. This does not limit the gener- 
ality of the conclusions obtained from such comparison due to the pronounced stress concentrations around 
the fiber/matrix interface, which give rise to the evolution of the local stress and inelastic strain fields whose 
complexity depends on the loading type. For instance, under axisymmetric loading, the inelastic strain field 
evolution depends only on the radial coordinate whereas under axial shear loading the angular dependence 
is also present, producing stress and inelastic strain fields that are two-dimensional. It is both the charac- 
ter and the magnitude of these local stress and inelastic strain fields for dilute fiber volume fractions that 
must be reproduced by the High-Fidelity Generalized Method of Cells, upon comparison with the multiple 
concentric cylinder model’s predictions, that will validate the theory’s predictive capability. Alternatively, 
at non-dilute fiber volume fractions, the presence of adjacent fibers will affect the local stress and inelas- 
tic strain fields. Therefore, the theory’s predictions are expected to reflect this interaction relative to the 
multiple concentric cylinder model which does not directly account for the adjacent fibers’ presence. Since 
the multiple concentric cylinder model cannot be used to generate a solution under transverse loading, the 
finite-element method was employed for this loading type to validate the theory’s predictive capability. 


Table 1 . Elastic and thermal parameters of the transversely isotropic graphite fiber. 


Ea (GPa) v a 

E t (GPa) 

v T G a (GPa) 

(10- 6 /°C) 

a T (10- 6 /°C) 

388.20 0.41 

7.60 

0.45 14.90 

-0.68 

9.74 


Table 2. Elastic, plastic and thermal parameters of the isotropic aluminum matrix. 


E (GPa) v 

G (GPa) 

a (10 - 6 /°C) G y (MPa) E s (GPa) 

72.40 0.33 

27.22 

22.5 286.67 11.70 


A graphite/aluminum (gr/al) system is considered with the thermal, elastic and plastic moduli of the in- 
dividual constituents given in Tables 1 and 2. The graphite fiber is linearly elastic and transversely isotropic, 
and the isotropic aluminum matrix is modeled using the classical incremental plasticity theory with isotropic 
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hardening. For this inelastic constitutive model, Mendelson’s method of successive iterations described 
briefly below, and in detail elsewhere (Mendelson, 1986; Williams and Pindera, 1997), is employed to solve 
the system of equations (45) for the unknown coefficients and the plastic strain distributions in 

the individual subcells. Towards this end, the plastic strain field in the (^ 7 ) subcell at the current load level 
is expressed in terms of the known initial distribution from the preceding loading state plus an increment 
that results from the imposed load increment, 

7) (sf + ,#) (56) 

The plastic strain increments at the individual locations within the subcell are calculated using the Prandl- 
Reuss flow rule expressed in terms of so-called modified total strain deviators e { - , rather than deviatoric 
stresses, as follows (omitting the superscript (^ 7 ) for notational clarity) 

/ 

c/ef,: = zr^de p (57) 

3 e eff 

where e V} = s l3 - 1/3 £kk^ij ~ < 7 previous, e e // = and the e ff ect i ve plastic strain increment 

de p is 

de p = e e ff — a/ (58) 

where /i is the shear modults. 

Without loss of generality in validating the theory’s predictive capability, we take the elastoplastic re- 
sponse of the aluminum matrix to be bilinear, with the effective stress < 7 ( 5 ?) given by 


a(s p ) = G y + H p e p 


(59) 


where a y is the yield stress in simple tension, and H p is the slope of the effective stress-plastic strain curve, 
related to the secondary modulus E s in the bilinear stress-strain representation of the elastoplastic response 
(see Table 2) as follows 


H p = 


EE S 

E-E s 


(60) 


where E is the Young’s modulus of the isotropic aluminum. 

The above form of the incremental plasticity equations is completely equivalent to the classical one, 
but has the added advantage of producing very quick convergence when used in conjunction with Mendel- 
son’s iterative scheme for the solution of Eq. (45), as demonstrated in our previous investigations. The 
implementation of these equations is facilitated by the following loading condition for plastic loading 


1 - 


a 

?>H e e ff 


>0 


(61) 


shown previously to produce very fast convergence. 

The fiber volume fractions ( Vf ) of the gr/al composite employed in the validation studies are 0.05 and 
0.25, respectively, which henceforth will be referred to as dilute and non-dilute fiber volume fractions. 
Figure 3 shows the volume discretizations of the repeating unit cell for these two cases. In both cases, the 
unit cells contain 36 x 36 subcells appropriately dimensioned to approximate the circular fiber shape with 
sufficient accuracy. This was accomplished using a commercial optimization code. 
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Figure 3. Volume discretization of the repeating unit cells employed in the analysis of unidirectional fiber- 
reinforced, metal-matrix composites: (left) Vf = 0.05; (right) Vf = 0.25. 


3.1 Axisymmetric thermal loading 

As the first step, we compare the predictions of the High-Fidelity Generalized Method of Cells and the 
multiple concentric cylinder model for the thermal response of the gr/al unidirectional composite during 
spatially uniform cooldown from 500° C to 25° C. This is the simplest case which results in axisymmetric 
displacement, stress and (plastic) strain fields within the individual phases of the multiple concentric cylinder 
model. As is well-known, the macroscopic axial response for this loading case can be predicted with good 
accuracy by the simplest micromechanics models irrespective of the fiber volume fraction. Thus a true test 
of the theory’s predictive capability is the evolution of the plastic strain field in the aluminum phase as a 
function of temperature, which should retain an axisymmetric character for the dilute (i.e., noninteracting) 
fiber volume fraction composite. This is the focus of this section. 

Figure 4 presents the macroscopic axial and transverse response of the gr/al unidirectional composite 
with the two fiber volume fractions as a function of temperature generated by the two models, which pro- 
vides information on the initiation of yielding in the aluminum phase. For the dilute fiber volume fraction 
composite, Fig. 4(top), yielding initiates slightly above 300° C, while in the non-dilute case this occurs at 
around 250° C, Fig. 4(bottom). In order to obtain initiation of yielding for the two cases at comparable 
temperature levels when cooled from a common temperature, the yield stress of the aluminum phase in 
the dilute fiber volume fraction composite was taken to be one fourth of the yield stress employed in the 
non-dilute fiber volume fraction composite reported in Table 2. 

As expected, both models predict virtually identical macroscopic thermal response in the axial direction 
for the dilute and non-dilute fiber volume fraction cases. The same observation holds true for the effective 
transverse response. Based on the macroscopic thermal response shown in Fig. 4(top), the effective plastic 
strain distributions for the dilute fiber volume fraction composite have been generated at 300°, 200° and 25° 
C using the High-Fidelity Generalized Method of Cells and the multiple concentric cylinder model. The 
effective plastic strain distributions are calculated by integrating the effective plastic strain increment de v at 
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each load increment along the entire loading history at the particular point, i.e., 

£ P = J d£ p 

where de p is given by Eq. (58). These distributions are shown in Fig. 5. In this dilute fiber volume 
fraction case, little interaction is expected between adjacent fibers in the periodic square fiber array model 
employed in the High-Fidelity Generalized Method of Cells. Thus the character of the effective plastic strain 
distribution in the aluminum phase predicted by this model should be essentially axisymmetric and therefore 
comparable to that predicted by the multiple concentric cylinder model. This is indeed the case with respect 
to both the distribution as well as the magnitude, as observed in Fig. 5. It is remarkable that the magnitude 
of the effective plastic strain field is predicted with high accuracy by the High-Fidelity Generalized Method 
of Cells even during the initial stages of yield initiation at 300° C shown in the top portion of Fig. 5. This 
demonstrates that the method is sufficiently sensitive to capture accurately the plastic strain field even when 
the effective plastic strain magnitudes are very small. 

In the case of the non-dilute fiber volume fraction composite, the effective plastic strain distributions 
have been generated at 200°, 100° and 25° C, based on the macroscopic thermal response given in Fig. 
4(bottom). These distributions are shown in Fig. 6. During the initial stages of yield initiation, the effective 
plastic strain distribution predicted by the High-Fidelity Generalized Method of Cells at 200° C retains an 
axisymmetric character and compares favorably with the multiple concentric cylinder model prediction. The 
interaction with adjacent fibers becomes apparent at the lower temperatures as the magnitude of the effective 
plastic strain predicted by the high-fidelity model within the plastic zone increases, losing its axisymmetric 
character away from the fiber. In the immediate vicinity of the fiber/matrix interface, however, the axisym- 
metric character of the effective plastic strain distribution is preserved and compares very well in magnitude 
with the multiple concentric model prediction. 
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Figure 4. Macroscopic thermal response of a unidirectional gr/al composite during spatially uniform cooldown 
from 500° C to 25° C: (top) Vf = 0.05; (bottom) Vf = 0.25. (Note that for Vf = 0.05, an altered matrix 
yield stress, a y = 286.67/4 MPa was employed, see Table 2). 
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Figure 5. Effective plastic strain distributions in the aluminum phase of a unidirectional gr/al composite with 
Vf = 0.05, cooled from 500° C, at 300° C (top), 200° C (middle), and 25° C (bottom). Comparison of the 
High-Fidelity Generalized Method of Cells predictions (left column) with the multiple concentric cylinder 
results (right column). 
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Figure 6. Effective plastic strain distributions in the aluminum phase of a unidirectional gr/al composite with 
Vf = 0.25, cooled from 500° C, at 200° C (top), 100° C (middle), and 25° C (bottom). Comparison of the 
High-Fidelity Generalized Method of Cells predictions (left column) with the multiple concentric cylinder 
results (right column). 
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3.2 Axial shear loading 

In comparing the predictions of both models for the axial shear loading case in the xi~xs plane (see Fig. 1), 
the following features of the local stress fields should be kept in mind. For the dilute fiber volume fraction 
composite in the elastic region, the cri3 shear stress distribution in the fiber region obtained from the multiple 
concentric cylinder model is uniform while the a 12 shear stress is zero. The evolution of plastic deformation 
in the matrix phase introduces slight nonuniformities in these distributions. This is also true in the non-dilute 
case since this model does not account directly for the interaction due to adjacent fibers. Further, the <713 
shear stress distribution is symmetric with respect to the horizontal and vertical planes through the fiber’s 
center, while the a 12 shear stress distribution is antisymmetric irrespective of the fiber volume fraction at all 
loading levels. With the exception of the interaction effect due to adjacent fibers in the non-dilute case, the 
above features of the local stress fields must be captured by the High-Fidelity Generalized Method of Cells. 

We start by presenting the macroscopic shear stress-strain predictions of the two models due to shearing 
by the applied average strain 613. These responses are given in Fig. 7 for both the dilute and non-dilute cases, 
and form the basis for the microscopic stress and effective plastic strain distributions presented at different 
applied shear strain levels in the subsequent figures. In generating the shear response, the properties of 
the aluminum matrix given in Table 2 were employed while the fiber was assumed to be isotropic with the 
shear modulus G = 137.7 GPa, following Williams and Pindera (1997). This produces a fiber/matrix shear 
modulus ratio of approximately 5, which provides a greater property mismatch than that obtained using the 
properties of the transversely isotropic graphite fiber given in Table 1 . 

In the case of the dilute fiber volume fraction composite, Fig. 7(top), the predictions of the High-Fidelity 
Generalized Method of Cells and the multiple concentric cylinder model are identical in both the elastic 
and elastoplastic region, with macroscopic yielding occuring around 613 = 0.3%. Based on the observed 
macroscopic response, the microscopic stress distributions are given at 613 = 0.1, 0.5, and 0.75%, while 
the microscopic effective plastic strain distributions are given at 613 = 0.3, 0.5, and 0.75%. The elastic 
shear stress distributions at 613 = 0.1% are included because the character of the microscopic stress field is 
affected by the growth of the plastic zone in the aluminum matrix even at low fiber volume fractions, unlike 
the axisymmetric thermal loading case discussed in the previous section. In the case of the non-dilute fiber 
volume fraction composite, Fig. 7(bottom), the predictions of the High-Fidelity Generalized Method of Cells 
and the multiple concentric cylinder model are also identical in the elastic region, with macroscopic yielding 
occuring around 613 = 0.225%. In the elastoplastic region, the macroscopic shear response predicted by 
the High-Fidelity Generalized Method of Cell is slightly lower relative to the multiple concentric cylinder 
model prediction, indicating that the interaction with adjacent fibers is becoming noticeable. Based on the 
observed macroscopic response, the microscopic stress distributions are given at 613 = 0.1, 0.5, and 0.75% 
as in the dilute case, while the microscopic effective plastic strain distributions are given at 613 = 0.25, 0.5, 
and 0.75%. 

The microscopic shear stress a is distributions predicted by the High-Fidelity Generalized Method of 
Cells and the multiple concentric cylinder model for the dilute case are compared in Fig. 8 at the given 
applied shear strain levels. Both the character and the magnitude of the local shear stress field is captured 
with very good accuracy by the high-fidelity model in this dilute case, including the initially nearly uniform 
shear stress distribution in the fiber and the highly nonuniform distribution in the matrix phase in the elastic 
region, and the effect of matrix plasticity on these distributions with progressive loading. Of particular 
significance is the symmetry of the a is field with respect to the horizontal and vertical planes through the 
fiber’s center during all stages of loading. The effect of plasticity on the character of the shear stress field 
is also dramatic, and perhaps more apparent, in the case of the a 12 distributions presented in Fig. 9 due 
to the antisymmetric character of the stress field with respect to the horizontal and vertical planes. In this 
case, yielding of the matrix phase skews the a 12 distributions away from the bisectors of the horizontal 
and vertical axes through the fiber’s center towards the vertical axis as observed in the middle and bottom 
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portions of Fig. 9. This phenomenon is very well captured by the High-Fidelity Generalized Method of 
Cells, as is the perfect antisymmetry of the a 12 field with respect to the horizontal and vertical planes. 
As in the preceding case, the magnitude of the a 12 field predicted by the high-fidelity model compares 
very well with the multiple concentric cylinder model prediction at the three applied shear strain levels. It is 
therefore not surprising that the effective plastic strain distributions predicted by the high-fidelity model also 
compare very well with the multiple concentric cylinder model predictions both in character and magnitude, 
as observed in Fig. 10. This is true even during the early stages of yield initiation when the effective plastic 
strain magnitudes are very small. 

In the non-dilute case, the effect of adjacent fibers on the microscopic shear stress fields generated by 
the High-Fidelity Generalized Method of Cells becomes apparent. The magnitude of this effect depends on 
the particular shear stress component and the extent of plastic deformation. In the case of the shear stress 
a 13 distributions shown in Fig. 11, the shear stress in the fiber is nearly uniform in the elastic region and 
compares very well in magnitude with the multiple concentric cylinder model result (Fig. 11 (top)). The 
matrix shear stress field is also well captured by the high-fidelity model. As the extent of the plastic zone 
grows, the interaction with adjacent fibers becomes more apparent both in the fiber and the matrix phases. 
Specifically, the high-fidelity model predicts greater extent of the fiber shear stress nonuniformity, and a 
smaller shear stress gradient in the matrix phase in the plane of shearing, relative to the multiple concentric 
cylinder results. 

In the case of the shear stress a 12 distributions shown in Fig. 12, the shear stress in the fiber becomes 
nonuniform already in the elastic region. However, the basic character and magnitude of the shear stress 
field is essentially the same as that predicted by the multiple concentric cylinder model. As in the dilute 
case, plastic deformation in the matrix phase skews the shear stress field away from the bisectors of the 
horizontal and vertical axes through the fiber’s center towards the vertical axis, albeit not to the same ex- 
tent as in the dilute case and as predicted by the multiple concentric cylinder model. In general, in the 
non-dilute case the overall character of the shear stress a 12 distributions in the elastoplastic region better 
resembles the corresponding multiple concentric cylinder model predictions than does the character of the 
a 13 distributions. 

The influence of adjacent fibers in the periodic array model employed by the High-Fidelity Generalized 
Method of Cells is also apparent in the effective plastic strain distributions shown in Fig. 13 for the non- 
dilute case. During the initial stages of yielding, this influence is not noticeable (Fig. 13 (top)), but increases 
at the higher levels of deformation (Fig. 13(middle and bottom)). Specifically, the constraint due to the 
adjacent fibers retards the growth of the plastic zone in the matrix phase in the direction normal to the 
plane of shearing, and enhances it in the shearing plane itself. In the immediate vicinity of the fiber/matrix 
interface, however, both the character and magnitude of the effective plastic strain field predicted by the 
high-fidelity model compare very well with the multiple concentric cylinder model prediction. 
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Figure 7. Macroscopic shear stress-strain response, a 13 vs 613 , of a unidirectional gr/al composite due to 
shearing in the x\ — X 3 plane: (top) Vf = 0.05; (bottom) Vf = 0.25. Comparison of the High-Fidelity 
Generalized Method of Cells predictions with the multiple concentric cylinder model results. (Note that in 
this case the graphite fiber is treated as isotropic with G = 137.7 GPa). 
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Figure 8. Shear stress a 13 distributions in the individual phases of a unidirectional gr/al composite with 
Vf = 0.05 at the applied average shear strain 613 of 0.1% (top), 0.5% (middle), and 0.75% (bottom). 
Comparison of the High-Fidelity Generalized Method of Cells predictions (left column) with the multiple 
concentric cylinder results (right column). Colorbar in MPa. (Note that the dot at the fiber’s center is a 
plotting artifact). 
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Figure 9. Shear stress a 12 distributions in the individual phases of a unidirectional gr/al composite with 
Vf = 0.05 at the applied average shear strain 613 of 0.1% (top), 0.5% (middle), and 0.75% (bottom). 
Comparison of the High-Fidelity Generalized Method of Cells predictions (left column) with the multiple 
concentric cylinder results (right column). Colorbar in MPa. 
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Figure 10. Effective plastic strain distributions in the aluminum phase of a unidirectional gr/al composite 
with Vf = 0.05 at the applied average shear strain 613 of 0.3% (top), 0.5% (middle), and 0.75% (bottom). 
Comparison of the High-Fidelity Generalized Method of Cells predictions (left column) with the multiple 
concentric cylinder results (right column). 
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Figure 11. Shear stress a 13 distributions in the individual phases of a unidirectional gr/al composite with 
Vf = 0.25 at the applied average shear strain 613 of 0.1% (top), 0.5% (middle), and 0.75% (bottom). 
Comparison of the High-Fidelity Generalized Method of Cells predictions (left column) with the multiple 
concentric cylinder results (right column). Colorbar in MPa. (Note that the dot at the fiber’s center is a 
plotting artifact). 
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Figure 12. Shear stress a 12 distributions in the individual phases of a unidirectional gr/al composite with 
Vf = 0.25 at the applied average shear strain 613 of 0.1% (top), 0.5% (middle), and 0.75% (bottom). 
Comparison of the High-Fidelity Generalized Method of Cells predictions (left column) with the multiple 
concentric cylinder results (right column). Colorbar in MPa. 
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Figure 13. Effective plastic strain distributions in the aluminum phase of a unidirectional gr/al composite 
with Vf = 0.25 at the applied average shear strain 613 of 0.3% (top), 0.5% (middle), and 0.75% (bottom). 
Comparison of the High-Fidelity Generalized Method of Cells predictions (left column) with the multiple 
concentric cylinder results (right column). 
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3.3 Transverse normal loading 

The commercial finite-element program ANSYS was employed to generate the response of the gr/al unidi- 
rectional composite under transverse loading. Refined representations of the repeating unit cell containing 
1,782 elements with 4,989 nodes for the dilute case, and 1,260 elements with 3,909 nodes for the non-dilute 
case, were constructed using 8-noded, plane-strain quadrilateral elements. External loading was specified 
in the form of uniform displacements of opposite signs but equal magnitudes applied to the horizontal and 
vertical boundaries of the unit cell under the condition of plane strain in the out-of-plane direction. This 
type of loading represents pure shearing in the coordinate system rotated by 45° about the fiber’s center, 
and provides a critical test of the High-Fidelity Generalized Method of Cells’ predictive capabilities under 
transverse loading. In the context of the finite-element analysis, it also eliminates the need for an iterative 
solution procedure when transverse loading is applied in just one direction in order to ensure that the average 
stress in the other direction remains zero. The same loading was employed to generate the High-Fidelity 
Generalized Method of Cells results (i.e., 622 = —633 with en = 0). 

For this loading case, the actual properties of the transversely isotropic graphite fiber were employed in 
the calculations. Therefore, the matrix is 3.58 times stiffer than the fiber in the loading plane (see Tables 1 
and 2), in contrast with the axial shearing case where the fiber was approximately 5 times stiffer than the 
matrix in the plane of shearing. The smaller transverse stiffness of the fiber relative to the matrix in the 
present case is expected to produce substantially different trends than those observed in the axial shearing 
case both at the macroscopic and microscopic scales, providing a completely different test for the High- 
Fidelity Generalized Method of Cells’ predictive capability. This is indeed the case as demonstrated first 
in Fig. 14 which provides the comparison between the high-fidelity model and finite-element predictions 
of the macroscopic stress-strain response for the two fiber volume fraction composites. Specifically, the 
stress-strain response for the non-dilute case is lower relative to the dilute case. Further, there is virtually no 
difference between the predictions of the two models for the two cases. 

Included in Fig. 14 are the corresponding predictions of the original Generalized Method of Cells for the 
two fiber fraction cases. As observed, the initial elastic response is very well captured by the Generalized 
Method of Cells despite the inherent absence of shear coupling in this simpler micromechanical model. The 
predictive capability of this model is also quite good in the elastoplastic region, but does depend on the fiber 
volume fraction as expected. The reason for the original Generalized Method of Cells’ ability to model the 
macroscopic inelastic response of metal matrix composites with good accuracy will become apparent upon 
examination of the internal stress fields. 

Based on the macroscopic stress-strain curves, local field quantities were generated at the applied 633 
strain levels of 0.4 and 1.0%. Figure 15 shows the effective plastic strain distributions generated by the 
High-Fidelity Generalized Method of Cells which can be compared with the corresponding distributions 
generated by the finite-element model given in Fig. 16. As observed, both the magnitudes and the character 
of the distributions are virtually the same at the two strain levels. In particular, for the dilute case, initiation 
of yielding at the fiber/matrix interface at locations coincident with the loading axes is correctly captured by 
the high-fidelity model, as is the spread of the plastic zone at the final macroscopic strain. Yielding in the 
non-dilute case initiates halfway between adjacent fibers along the line rotated 45° about the fiber’s center, 
in contrast with the dilute case. This too is captured by the high-fidelity model with excellent accuracy, as 
is the effective plastic strain distribution at the final macroscopic strain. 

In order to highlight the differences between the High-Fidelity Generalized Method of Cells and the 
original Generalized Method of Cells, Fig. 17 illustrates the inplane shear stress a 23 distributions generated 
using the high-fidelity model at the same applied 633 strain levels as those presented in Fig 16. Clearly, the 
magnitude of the inplane shear stress is significant and cannot be neglected. The corresponding distributions 
obtained from the finite-element analysis are shown in Fig. 18. Comparison of the two sets of predictions 
reaffirms the predictive capability of the high-fidelity model. In stark contrast, the original Generalized 
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Method of Cells predicts that the inplane shear stress distributions are identically zero, which is a direct 
consequence of the absence of shear coupling. Apparently, the absence of shear coupling, while resulting in 
poor prediction of local distributions of certain stress components, has little effect on the method’s ability to 
accurately capture the macroscopic behavior. 

In order to explain the original Generalized Method of Cells’ capability to capture the macroscopic 
response of metal matrix composites in the inelastic region with good accuracy, it is neccessary to examine 
those stress distributions that are directly responsible for the inelastic behavior. In the case of the classical 
incremental plasticity theory employed herein, it is necessary to examine the effective stress distributions. 
These are illustrated in Figs. 19, 20, and 21 for the two fiber fraction cases generated at the same applied 
macroscopic strain levels as above using the High-Fidelity Generalized Method of Cells, finite-element 
analysis, and the original Generalized Method of Cells, respectively. Comparison of Figs. 19 and 20 leads 
to the same conclusions regarding the predictive capability of the high-fidelity model as the comparison 
of the effective plastic strain distributions shown in Figs. 15 and 16. Examination of the effective stress 
distributions generated using the original Generalized Method of Cells, Fig. 21, leads to the conclusion that 
the essential features of these distributions are generally the same as those obtained from the high-fidelity 
model and the finite-element analysis, Figs. 19 and 20, even though the local details differ somewhat. This 
explains the original Generalized method of Cells’ ability to model the elastoplastic response of metal matrix 
composites with sufficient accuracy. 

Finally, Figs. 22, 23, and 24 present the corresponding hydrostatic stress distributions generated using 
the High-Fidelity Generalized Method of Cells, finite-element analysis, and the original Generalized Method 
of Cells, respectively. As in the preceding cases, comparison of Figs. 22 and 23 indicates that the high- 
fidelity model captures both the magnitude and local distributions of the hydrostatic stress with very good 
accuracy relative to the finite-element predictions. In contrast, the accuracy with which these distributions 
are captured by the original Generalized Method of Cells is quite poor. It is only in a gross average sense 
that these distributions are accurate over large areas of the repeating unit cell, precluding the possibility of 
accurately modeling local matrix degradation or failure due to the hydrostatic stress component. 
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Figure 14. Macroscopic transverse stress-strain response, <733 vs 633, of a unidirectional gr/al composite 
due to loading in the X 2 — £3 plane by 622 = —633 with en = 0: (top) Vf = 0.05; (bottom) Vf = 0.25. 
Comparison of the High-Fidelity Generalized Method of Cells predictions with the results from the finite- 
element analysis and the original General Method of Cells. 
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Figure 15. Effective plastic strain distributions in the aluminum phase of a unidirectional gr/al composite at 
the applied average transverse strain 633 of 0.4% (left column) and 1.0% (right column) obtained from the 
High-Fidelity Generalized Method of Cells: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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Figure 16. Effective plastic strain distributions in the aluminum phase a unidirectional gr/al composite at 
the applied average transverse strain 633 of 0.4% (left column) and 1.0% (right column) obtained from the 
finite-element analysis: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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Figure 17. Inplane shear stress a 23 distributions in the individual phases a unidirectional gr/al composite at 
the applied average transverse strain €33 of 0.4% (left column) and 1.0% (right column) obtained from the 
High-Fidelity Generalized Method of Cells: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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Figure 18. Inplane shear stress a 23 distributions in the individual phases a unidirectional gr/al composite at 
the applied average transverse strain 633 of 0.4% (left column) and 1.0% (right column) obtained from the 
finite-element analysis: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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Figure 19. Effective stress distributions in the individual phases a unidirectional gr/al composite at the 
applied average transverse strain 633 of 0.4% (left column) and 1.0% (right column) obtained from the 
High-Fidelity Generalized Method of Cells: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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Figure 20. Effective stress distributions in the individual phases a unidirectional gr/al composite at the 
applied average transverse strain 633 of 0.4% (left column) and 1.0% (right column) obtained from the 
finite-element analysis: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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Figure 21. Effective stress distributions in the individual phases a unidirectional gr/al composite at the 
applied average transverse strain 633 of 0.4% (left column) and 1.0% (right column) obtained from the 
original Generalized Method of Cells: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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Figure 22. Hydrostatic stress distributions in the individual phases a unidirectional gr/al composite at the 
applied average transverse strain 633 of 0.4% (left column) and 1.0% (right column) obtained from the 
High-Fidelity Generalized Method of Cells: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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Figure 23. Hydrostatic stress distributions in the individual phases a unidirectional gr/al composite at the 
applied average transverse strain 633 of 0.4% (left column) and 1.0% (right column) obtained from the 
finite-element analysis: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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Figure 24. Hydrostatic stress distributions in the individual phases a unidirectional gr/al composite at the 
applied average transverse strain 633 of 0.4% (left column) and 1.0% (right column) obtained from the 
original Generalized Method of Cells: (top) Vf = 0.05; (bottom) Vf = 0.25. 
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4 Summary and Conclusions 


A recently developed theory for periodic multiphase materials, called High-Fidelity Generalized Method 
of Cells, which was previously limited to thermoelastic phases has been extended herein to admit inelastic 
constitutive behavior of the individual phases. This extension has been carried out in a general fashion in 
order to enable the use of different constitutive models for the phase behavior, including classical incre- 
mental plasticity theory, creep models, as well as various unified viscoplasticity theories. The theoretical 
framework combines elements of the homogenization technique with the higher-order theory for function- 
ally graded materials developed previously by the authors. These features provide a basis for consistent 
approximation of the displacement field at the local level together with consistent application of periodic 
boundary conditions imposed on the deformation of the repeating unit cell, which characterizes the mater- 
ial’s microstructure. Further, the higher-order displacement field approximation at the local level employed 
in the present approach provides the necessary coupling between the local normal and inplane shear stress 
fields and the macroscopically applied loading. This coupling dramatically improves the accuracy of esti- 
mating the local stress and inelastic strain fields relative to the original Generalized Method of Cells which 
is based on a first-order displacement approximation at the local level. The high-fidelity model’s predictive 
capability to capture both the macroscopic response and the local stress and inelastic strain fields has been 
demonstrated through comparison with the results of two analytical solutions and finite-element element 
analysis of the inelastic response of a unidirectional gr/al composite, based on the incremental plasticity 
theory, subjected to different types of loading. 

The primary result of the High-Fidelity Generalized Method of Cells is the closed- form constitutive 
equation for the macroscopic thermoinelastic response of multiphase materials, subjected to arbitrary mul- 
tiaxial macroscopic thermomechanical loading, which possess microstructures characterized by repeating 
unit cells with arbitrary reinforcement distributions. This is a direct result of the use of periodic bound- 
ary conditions that follow from the homogenization approach’s framework. Thus the macroscopic inelastic 
response can be generated irrespective of whether or not a repeating unit cell possesses planes of mater- 
ial symmetry. This can be done easily for any combination of macroscopically applied thermomechanical 
loads, in contrast with the standard finite-element analyses of periodic composites. The repeating unit cell’s 
construction is simple due to the employed volume discretization that produces a rectangular grid whose 
subcells are appropriately assigned different material properties and dimensions so as to mimic a multiphase 
material’s actual microstructure. Further, the computational speed with which the macroscopic inelastic 
response and local stress and plastic strain fields are generated is sufficiently fast (typically on the order of 
a few minutes on the DEC-Alpha DS20E 6/667 machine for the investigated cases) for reasonably detailed 
volume discretizations of a repeating unit cell. These features of the presented theory facilitate investiga- 
tions of the impact of different materials architectures on both the macroscopic and local responses, taking 
into account inelastic phase behavior, in an efficient and accurate manner. They also make it straightforward 
to incorporate the theory into a structural analysis computer code as a subroutine. 
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Appendix: Multiple Concentric Cylinder Model 

A brief outline of the solutions for the displacement field, from which the corresponding stress field can 
be generated, in a unidirectional elastoplastic composite subjected to axisymmetric and axial shear loading 
is given in this section based on the multiple concentric cylinder model. The solutions are fully analytical 
and therefore provide a good basis for comparison with the predictions of the outlined higher-order theory 
despite the differences in the model geometries. The cylindrical coordinate system x—r—6 (where x denotes 
the direction along the cylinder’s axis) is used to formulate the problem and solve for the displacement and 
stress fields in a multilayered cylinder consisting of an elastic core and an arbitrary number of fully bonded 
inelastic concentric shells. The elastic core is denoted by the superscript 1 and the outermost cylindrical 
shell by n. The inner radius of the fcth shell is denoted by rk-i and the outer radius by r&. 

Axisymmetric loading 

Under axisymmetric loading that may involve a combination of a macroscopically uniform axial stress 
or strain, uniform temperature field and biaxial tension/compression, the displacement field in the individual 
layers of a multiple concentric cylinder has the following form 

u x = s° xx x, u r = u r (r), u e = 0 (Al) 


where £ xx is the uniform axial strain component. Therefore, the equilibrium equations in the individual 
layers expressed in terms of displacements reduce to the single ordinary differential equation 


d 2 U r 1 du r _ _ 1 (Cri — Cpi) ins \ I Cri d£™(r) 

dr 2 r dr r 2 r ^ C rr £n V ^ C rr dr 

i=x,Q,r i=x,Q,r 


(A2) 


where s™{r) are the inelastic strain distributions that depend implicitly on the radial displacement field 
u r (r), and C%j (i,j = x,0,r) terms are the components of the elastic stiffness matrix expressed in the 
cylindrical coordinate system. For the elastic core, the right hand side of the above equation becomes 
zero. In the present situation, we consider thermal loading due to a spatially uniform temperature change. 
Therefore, the solution to the above equation in each layer of the composite cylinder is obtained subject to 
the boundary conditions 

Kr( r n) = 0, (A3) 

the interfacial displacement and traction continuity conditions 

v%- 1 {r k -i) = u*(r k -i), a k r ~ l {r k -i) = a k r (r k -i) (A4) 


where k 
blage 


= 2 and the axial equilibrium condition for the entire multiple concentric cylinder assem- 



(A5) 


where A c is the cross-sectional area of the assemblage. 

Using standard techniques, the solution of the equilibrium equation in each layer is obtained in the form 


u r (r) 


, A 2 1 
Ai r + — + — 
r 2r 


E 


( Cm + Cei ) 


C r 


£ u( r ')r' dr + 


r k-l l— x Q r 


r k ~ 1 i=x,6,r 

(C r i — C$i) j , dr 
r 


E \ E 7£4?<r»-,)r<^-l) 


i=x,0,r 


C r 

C~ r 


(A6) 
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where rk-i < r < r&. The unknown coefficients Ai,A 2 , axial strain £® x , and inelastic strain distributions 
in the core and the fcth shell are obtained at each increment along the imposed thermal loading path by 
following the solution procedure outlined by Pindera et al. (1993). The imposed spatially uniform thermal 
loading appears in the solution procedure through the application of the interfacial traction continuity and 
external boundary conditions upon expressing tractions in terms of strains (and thus the displacement field) 
using Hooke’s law. 

Axial shear loading 

Under axial shear loading by homogeneous displacements or tractions that produce a uniform axial shear 
strain in an equivalent homogenized medium, the displacement field takes the following form, referred to 
the coordinate system of Fig. Al, 

u x {r,6) = 4>(r,0) — streets 9, u r (x,9) = £i 2 xcos9 1 uq{x,9) = —e^xsinO (A7) 


where e\ 2 is the uniform shear strain in an equivalent homogenized medium generated by homogeneous 
displacement or traction boundary conditions in the horizontal plane. The function </>(r, 9) represents the 
deviation in the axial shear deformation of a layer from the solution for the equivalent homogenized material 
with homogeneous properties and, in the presence of inelastic effects, depends implicitly on the inelastic 
strain field. Therefore, the equilibrium equations in the individual layers expressed in terms of displacements 
result in the following partial differential equation for the unknown function </>(r, 9) in each layer 


1 d dej) 1 d 1 2 (j) 
r dr T dr + r 2 d9 2 


l d(re™) 1 de% 
r dr r d9 


e(r, 9) 


(A8) 


The solution to the above equation is obtained subject to the homogeneous displacement boundary condi- 
tions imposed on the surface of the nth shell 


u x (r, 9) = £ 12 bcos9, u r {x,0) = £i 2 xcos9, uq{x,9) §= —£ 12 xsin.9 (A9) 


where b is the outer radius of the entire cylindrical assemblage, and the interfacial displacement and traction 
continuity conditions. For the above displacement field representation, these reduce to 

Ux( r k-i,0), (?xr \r k - 1 ,9) m\cr k xr {r k -i : 9) (A10) 


since the continuity of the radial and tangential displacement components u r (x, 9) and uq (x, 9) is identically 
satisfied, and the traction components a rr and cr r Q vanish. 

The solution to the partial differential equation governing the function </>(r, 9) is obtained in the form 


<f>{r, 0) = 


1 1 X °° A 
—=R 0 (r) + —=y^Rn(r) cos n9 

V27 t v 71 ” i 

v v n = 1 


where 


for n — 0, and 


Ro(r) = [Fio + f sf 0 (s)ds] Inr + F 20 ~ j shisfo(s) 
Jr k _ 1 Jr k - 1 


R n (r) = [F ln + T f s n+1 f n (s) ds]r n + [F 2n - T f s n+1 f n (s) ds]r~ 
J r k — 1 J r k — 1 


(All) 


(A 12) 


(A 13) 


for n > 0. The functions / n (r) are the coefficients of the Fourier series representation of e(r, 9). They are 
given by 

/ n (r) = —j= f e(r, 9) cosn9d9 (A14) 


7 r 
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The unknown coefficients F\ n and F 2 n appearing in the solution for 0) in each layer are obtained by 
applying the interfacial displacement and traction continuity conditions, external boundary conditions, and 
employing an iterative solution procedure at each load increment in the manner described by Williams and 
Pindera (1997). Convergent solutions in the presence of inelastic effects are typically obtained with 25 — 30 
harmonics in the Fourier series representation of </>(r, 6). 


EQUIVALENT HOMOGENIZED MEDIUM 


Displacement bc's for 
homogeneous shearing 
in the x,-x 3 plane 

^1 “ ^ 12^2 

u 2 £ 12^1 

u 3 = 0 


Corresponding 
displacement bc's 
on the lateral surface 
of a homogenized 
cylinder 

u x = e° 2 bcos0 
u r = e° 2 xcos0 
u„= e° 2 xsin0 



x„x 


Figure Al. Boundary conditions for the multiple concentric cylinder model subjected to axial shearing in 
the x\ — X 2 plane. 
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